Sliding elastic lattice: an explanation of the motion of superconducting vortices 
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We introduce a system where an elastic lattice of particles is moved slowly at a constant velocity 
under the influence of a local external potential, construct a rigid-body model through simplification 
processes, and show that the two systems produce similar results. Then, we apply our model to 
a superconducting vortex system and produce path patterns similar to the ones reported in [Lee 
et al., Phys. Rev. B 84, 060515 (2011)] suggesting that the reasoning of the simplification processes 
in this paper can be a possible explanation of the experimentally observed phenomenon. 



The systems having many interacting particles, such 
as superconducting vortices [1], charge density waves 
[2], colloids [3], and Wigner crystals [4], have been of 
great interest in condensed matter physics because they 
show rich dynamic and static properties and the various 
phases that originate from the noncquilibrium dynam- 
ics and statistics. For example, the motion of super- 
conducting vortices under the influence of pinning po- 
tentials has been extensively studied in experiments [5- 
11], theories [12, 13], and simulations [14, 15], because 
the vortex system can be a model for the study of the 
dynamics of many interacting particle systems. More- 
over, understanding the motion of vortices is important 
for the development of the superconductivity-based de- 
vices [16]. One of the main topics of the vortex dynamics 
study is how to effectively suppress the motion of vor- 
tices because the heat dissipation caused by the motion 
of vortices can destroy the superconducting state. There- 
fore, extensive studies [5-15] have focused on the pinning 
effects on vortex motions when the vortex-pinning inter- 
actions are large enough to disturb the regularity of the 
Abrikosov vortex lattices [1]. However, vortex dynamics 
under weak pinning potentials has been rarely studied to 
date. Here, we ask ourselves, how do the moving vortices 
behave when the vortex- vortex interactions dominate the 
vortex-pinning interactions? 

Recently, the motion of the superconducting vortices 
in a single crystal NbSe2 sample was observed over 7 days 
using a scanning tunneling microscopy when the vortices 
were driven at a few pm/s by a slowly decaying mag- 
netic field at the rate of a few nT/s [17]. The paper 
emphasizes two peculiar features in the observed vortex 
dynamics: the nonuniform motion of the vortices and the 
existence of favored locations by the moving vortices. A 
2-dimensional molecular dynamics simulation is a good 
candidate for the understanding of the observed features. 
However, it requires a long computation time with a large 
number of particles. Indeed, an extensive 2-dimcnsional 
molecular dynamics simulation with about 8000 particles 
has been done [18], but the simulation was not able to 
provide a physical insight into the observed peculiar fea- 
tures due to the limitation of the computing power. In 



this paper, we propose a model that can provide a good 
physical insight into the origin of the observed features 
without use of extensive computing power. 

The first goal of this paper is to study the properties 
of a 2-dimensional elastic lattice and construct a compu- 
tationally efficient approximate model that captures the 
basic mechanism of the behavior of the slowly moving 
clastic lattice under the influence of a weak local attrac- 
tive potential fixed in space. The second goal of the paper 
is to show that this model can be applied to a supercon- 
ducting vortex system under certain conditions and pro- 
duce similar result, thus providing a possible explanation 
of the experimentally observed features of superconduct- 
ing vortices reported in [17]. 

A superconducting sample develops current vortices 
when an external magnetic field is applied to the sample 
perpendicular to the surface of the sample. The vortices 
exert repulsive forces among them due to electromagnetic 
forces. Vortices can be treated as classical particles with 
repulsive forces among them packed inside the boundary 
of the sample. They form a triangular lattice structure 
when the external magnetic field is strong enough. As the 
strength of the magnetic field decreases, the vortices near 
the boundary of the sample disappear and the number 
of vortices in the sample decreases making the distances 
among them increase. This causes vortex motion roughly 
in radial directions away from the center of the sample. 
It was reported that the vortices move non-uniformly at 
the speed on the order of pm/s with the triangular lattice 
structure maintained in the scanning window and there 
are favored locations where vortices visit and stay longer 
than at other locations [17]. 

To make a model that produces the observed behav- 
iors, we assume a quasistatic equilibrium because the vor- 
tices move extremely slowly. We also assume there are 
impurities in the sample introducing external attractive 
potentials and that they are weak because of the pristine 
nature of the sample. We are interested in the motion 
of the vortices near the impurity close to the scanning 
area. To estimate the computational cost of the direct 
simulation, consider particles between 2 fixed boundary 
points at the origin and I in a 1-dimensional space and 
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the motion of the particles as the particles next to the 
right boundary are removed one by one. Removing 1 par- 
ticle causes a particle at x move about x/L of the lattice 
constant. The observation location needs to be less than 
L/10 so that the paths of the particles at the observation 
location have at least 10 steps as the particles move for 
about 1 lattice constant. We need to remove at least 10 
particles so that the particles at the observation location 
move for at least 1 lattice constant. We need to have at 
least 100 initial particles so that the lattice constant does 
not change much as reported in [17] as 10 particles are 
being removed. In a 2-dimensional case, we need to have 
at least 100 x 100 particles making the direct simulation 
computationally too expensive to perform. 

We note that the direct cause of the motion of the vor- 
tices is the moving neighboring vortices. If we consider a 
collection of vortices of a certain area, then the motion of 
those vortices is caused by the repulsive forces from the 
surrounding moving vortices. Replacing the surrounding 
vortices with a repulsive boundary moving at a constant 
velocity to approximate the driving mechanism and in- 
troducing a weak external attractive potential fixed in 
space for the effect of the weak pinning caused by an 
impurity leads to a moving boundary model. However, 
this model did not produce the desired features possibly 
because we could not simulate enough particles limited 
by computing power, assuming the reported features in 
[17] are the results of a large number of vortices. This as- 
sumption is explained later in the simplification process 
from the elastic-lattice model to the rigid-body model. 

Assuming the pinning is weak enough not to break the 
lattice structure of the vortices, we replace the particles 
with repulsive forces in the previous model with a lat- 
tice of particles and springs. This motivates the sliding 
clastic lattice model. Consider a 2-dimensional triangu- 
lar lattice structure of particles with each connected to 6 
neighbors by springs of natural length 1 inside a certain 
boundary under the influence of a weak local external 
potential at a fixed location. We drag the boundary at 
a constant velocity slowly enough to ignore the masses 
of the particles and use quasistatic equilibrium analysis. 
We assume that the external potential is weak enough 
not to distort the lattice structure too much and that 
the external potential is short-ranged compared to the 
unit length of the lattice. We are interested in the paths 
of the particles near the external potential. 

Even if the elastic lattice model is computationally 
more efficient than the moving boundary model, the com- 
putational cost of the elastic lattice model is still high. 
To make a computationally more efficient model, we draw 
a virtual boundary containing several particles of our in- 
terest near the external potential resulting in the system 
divided into the area of our interest and the bigger envi- 
ronment containing the rest of the particles. There are 
three orders of strength involved in this picture. Firstly, 
the effective spring constant of the environment is much 



less than the spring constant of the individual springs be- 
cause the environment contains many springs connected 
in part in series and the combined spring constant de- 
creases when springs are connected in series. Secondly, 
the strength of the external potential is weak enough not 
to distort the lattice structure much, if we assume the 
external potential is weak compared to the strengths of 
the individual springs. Thirdly, the external potential is 
strong enough to disturb the locations of the particles 
of our interest because the effective spring of the envi- 
ronment is soft. Thus, we expect that physics does not 
change much when we replace the springs and the par- 
ticles inside the virtual boundary with a rigid body of 
those particles and replace the environment with a soft 
representative spring connected to the rigid body. 

We need to characterize the properties of the lattice 
before we simplify the system. To do that, we consider 
an elastic lattice with one particle located at the origin 
and the triangular structure oriented so that one of the 
6 directions of the edges is parallel to the x-axis, ignor- 
ing the existence of the external potential at the mo- 
ment. We cut out the lattice outside the circle of radius 
r centered at the origin when the lattice is at equilibrium 
with lattice constant 1. We call the particles that have 
all 6 neighboring particles inside particles, and the parti- 
cles that have less than 6 neighboring particles boundary 
particles. We fix the relative positions of the boundary 
particles to each other, letting the inside particles free to 
move. 

A 2-dimensional lattice is somewhat different from a 
1-dimcnsional lattice, because the 2-dimensional lattice 
becomes unevenly distorted when the center particle is 
moved by an external force while the 1-dimcnsional lat- 
tice becomes evenly compressed or stretched resulting in 
the effective spring constant of the system being inversely 
proportional to the size of the system. Thus, we need to 
characterize the distortion of the lattice. To characterize 
the distortion, we moved the center particle of the lattice 
to the right by 0.5 distance keeping the boundary parti- 
cles at the same positions. We found that the resulting 
displacements of the particles are quite local near the 
center, and the distortion of the hexagon of the 6 neigh- 
boring particles surrounding the center one is not negligi- 
ble. Thus, replacing the lattice near the center one with 
a rigid body is not reasonable in this situation. But if the 
functional form of the potential energy of each individ- 
ual spring has a shape that introduces more penalty away 
from the equilibrium than a quadratic function (obeying 
Hooke's law) does, then we can expect that the displace- 
ments become more spread and the distortion near the 
center decreases. To quantify the distortion, we define 
distortion D n of the hexagon surrounding the n-th par- 
ticle as 

n 2 _ St=l [ i x n,i — x' n i ) 2 + (y n ,i — y'n^j) 2 ] . . 

n ~ 6 ' [> 
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FIG. 1. (color online) (a) Distribution of distortions. The ra- 
dius of the boundary is 100. The values with quartic potential 
springs were multiplied by 6.22. (b) Maximum distortion vs. 
lattice size (radius). The vertical axis is the maximum dis- 
tortion occurring at the center particle moved by the external 
force, (c) force vs. lattice size (radius). The vertical axis is 
the force necessary to move the center particle by 0.5 distance. 
The forces with quartic potential springs were multiplied by 
4.95. (d) Force vs. Ax. The quartic potential springs are 
used. The radius of the lattice is 100. The hollow blue dots 
are the measured forces necessary to move the center particle 
by Ax. The solid curve is a fit by cubic function 5.72 x 10 _5 a; 3 . 



where x n ^ and y n ,i are the x and y positions of the ith 
neighboring particle of the distorted hexagon and x' n i 
and y' n i are the x and y positions of the ith neighboring 
particle of an undistorted hexagon. Figure 1(a) shows the 
distortions as a function of distance from the center with 
springs of V = ^(Ax) 2 potential energy (quadratic) and 
V = j{Ax) A potential energy (quartic) when the radius 
of the boundary is 100. The figure shows that the dis- 
tortions are more spread throughout the lattice and the 
maximum distortion occurring at the center is smaller 
with quartic potential springs than with quadratic poten- 
tial springs. Figure 1(b) shows the maximum distortion, 
which occurs at the center particle that was moved by 
the external force, as a function of the lattice size (the 
radius of the boundary) with quadratic potential springs 
and with quartic potential springs. The figure shows that 
the maximum distortion decreases faster with quartic po- 
tential springs than with quadratic potential springs as 
the lattice size increases implying the rigid body model 
is a better approximation with quartic potential springs 
than with quadratic potential springs given the same lat- 
tice size. 

Next, we characterize the effective spring of the lattice. 
Figure 1 (c) shows the amount of external force necessary 
to move the center particle by 0.5 distance as a function of 
the lattice size with quadratic potential springs and quar- 
tic potential springs. The figure shows that the force nec- 
essary to move the center particle decreases much faster 



with quartic potential springs than with quadratic poten- 
tial springs as the lattice size increases, and equivalently 
that the effective spring the lattice becomes softer much 
faster with quartic potential springs than with quadratic 
potential springs as the lattice size increases. 

So far, we compared the properties of the lattice of 
quadratic springs and the lattice of quartic springs. We 
also found that the results extend further when the higher 
exponents than 2 or 4 are used in the potential function 
of the springs. For the use in the rigid body model, we 
measured the spring constant k of the lattice of quartic 
springs when the radius of the lattice is 100 by measuring 
the forces necessary to move the center particle by Ax 
ranging from to 0.5. Figure 1(d) shows the force as 
a function of Aar. We find the spring constant k of the 
lattice to be 5.72 x 10~ 5 from the figure. 

Now, we compare the elastic lattice model and the 
rigid body model. The elastic lattice model is the same 
as described before with a lattice of the radius of 100, 
springs of V = j(Ax) A , and an external potential of 2- 
dimensional inverted Gaussian with standard deviation 
w/%/2 of 0.2/V2 and height h of 1.144 x 10~ 7 fixed at 
the origin. We rotated the lattice by arctan(\/3/5) ra- 
dian counterclockwise to produce a nontrivial path pat- 
tern, moved the boundary particles to the right in small 
intervals keeping the shape of the boundary the same 
but allowing the inside particles move freely, and then 
numerically calculated the equilibrium positions of the 
particles at each step by finding the positions of the in- 
side particles that minimize the total potential energy 
defined as 
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where x and y are vectors with each component rep- 
resenting x and y positions of each particle, N s is the 
number of springs of the system, and l n is the length 
of the nth spring that is determined by the positions of 
the two particles on the ends, and iV is the number of 
the inside particles. Figure 2(a) shows the path pattern 
near the external potential. The small black dots are 
the positions of the particles at each traveled distance by 
the boundary particles. The solid curves connecting the 
small black dots are the paths of each particle. The big 
gray dots are the positions of the particles with the center 
particle at the origin when there is no external potential. 
The figure shows that the particles tend to stay longer 
near the gray areas. We did a similar simulation with a 
rigid body of a rotated hexagonal shape of radius 20 and 
an effective spring of spring constant k of 5.72 x 10~ 5 as 
determined in Fig. 1(d) attached to the center particle. 
We moved the equilibrium position of the effective spring 
to the right in small intervals, and calculated the position 
of the center particle that minimizes the total potential 





FIG. 2. Paths of the particles near the external potential of 
(a) the elastic lattice model (b) the rigid body model when the 
system moves to the right in small intervals with the external 
potential located at the origin. The small black dots are the 
positions of the particles at each distance step, and the big 
gray dots are the positions of the particles with the center 
particle at the origin and no external potential. 



energy defined as 
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V( X ,y) = ^[( X -Xr + y2] 2 

// — .I 

(3) 

where x and y are the horizontal and vertical positions 
of the center particle of the rigid body, X is the hori- 
zontal equilibrium position of the center particle due to 
the effective spring, N is the number of particles of the 
rigid body, and x n and y n are the fixed horizontal and 
vertical positions of the nth particle in the rigid hexagon 
with respect to the center particle. The second term of 
the potential can be viewed as the potential of the center 
particle due to multiple Gaussian potentials, so the sys- 
tem can be viewed as one particle attached to a spring 
and moving through multiple external potentials located 
at the gray areas in the figure. Figure 2(b) shows the 
result. Figs. 2(a) and 2(b) show similar path patterns 
implying that the rigid body model is a good approxima- 
tion of the elastic lattice model. 

To produce the path pattern of the superconducting 
vortices reported in [17], we used a rigid hexagon of ra- 
dius 20 rotated by -arctan(v / 3/ll) radian, a spring of 
spring constant 1 attached to the center particle, and an 
external potential defined as 



G(x,y) 



■ _ Ae -{v?/*l + W/vl)^ for u > o 
_ Ae -(u 2 /<yl_+v 2 /al) otherwise 



(4) 



u = x cos ft — y sin ft 
v = xsinO + y cos ft, 

where A = 0.03, <r u+ — 0.6, a u - — <r v = 0.2, and 
ft = —132°. We moved the equilibrium position of the 
spring to the left. Figure 3 shows the path of the center 
particle and the overlaid paths of several particles around 
the center. This figure and Fig. 2 of [17] show similar 




FIG. 3. The path patterns of the rigid body model, (a) The 
path of the center particle, (b) The overlaid paths of several 
particles around the center particle. The lattice was rotated 
by -arctan(y / 3/ll) radian. An asymmetric Gaussian external 
potential was used. The equilibrium position of the effective 
spring was moved to the left in small intervals. 



path patterns even if the paths from our model make 
sharper changes of directions than the paths from the ex- 
perimental observation do. From this result, we propose 
that the pinning effects of the impurities in the sample 
used in the experiment of [17] were weak enough not to 
destroy the lattice structure of the vortices around them, 
but the large number of vortices in the sample made the 
effective spring of the vortex system soft enough to pro- 
duce the observed vortex dynamics. 

In summary, we went through several simplification 
processes of the superconducting vortex system introduc- 
ing the sliding elastic lattice model with a weak exter- 
nal potential and a rigid body model that approximates 
the behavior of the particles near the external potential, 
and showed similar path patterns from the elastic lattice 
model and the rigid body model. Then, we applied the 
rigid body model to the superconducting vortex system 
under a slowly decaying magnetic field, and showed that 
the path pattern from our model is similar to the ex- 
perimentally observed one suggesting that the reasoning 
leading to our model can be a possible explanation of the 
observed vortex dynamics reported in [17]. 
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